home *** CD-ROM | disk | FTP | other *** search
/ Practical Algorithms for Image Analysis / Practical Algorithms for Image Analysis.iso / CH_5.6 / fitspline / SPLINE.C < prev    next >
Encoding:
C/C++ Source or Header  |  1999-09-11  |  3.4 KB  |  96 lines
  1. /* SPLINE:      function fits splines to data points 
  2.  *                    usage: spline (nodes, iNode, image, width, height,
  3.  *                                                      gran, splineType);
  4.  *              The cardinal spline (splineType = 2) is an interpolating
  5.  *              spline that goes through all points. The B-spline is an
  6.  *              approximating spline that may not go through any points,
  7.  *              but gives a smoother fit.               
  8.  *              NOTE that the user must send this function a list of
  9.  *              points where the first and last points have been 
  10.  *              replicated (1st pt, 1st pt, subsequent pts, last pt, last pt).
  11.  */
  12.  
  13. #include <images.h>
  14.  
  15. #define INIMAGE(PT1,PT2) (PT1.x >= 0 && PT1.x < imgSize.x \
  16.               && PT1.y >= 0 && PT1.y < imgSize.y \
  17.               && PT2.x >= 0 && PT2.x < imgSize.x \
  18.               && PT2.y >= 0 && PT2.y < imgSize.y)
  19.  
  20.  
  21. long drawline8 (unsigned char **, struct point, struct point,
  22.                 struct point, unsigned char);
  23. long
  24. spline (nodes, n, image, imgSize, nT, splineType)
  25.      struct point *nodes;       /* list of node coord.s for spline */
  26.      long n;                    /* total number of node points */
  27.      unsigned char **image;     /* image in which to draw spline */
  28.      struct point imgSize;      /* size of image */
  29.      long nT;                   /* granularity of points between nodes */
  30.      short splineType;          /* =1 for approx spline; or =2 for interp */
  31.  
  32. {
  33.   struct point pt1, pt2;
  34.   long nM2;                     /* number of input pts minus 2 */
  35.   double deltaT;
  36.   double div;                   /* divisor */
  37.   long i, j;
  38.   double t, t1, t2, t3, t4;
  39.   double xOut, yOut;
  40.  
  41.   deltaT = 1.0 / nT;
  42.   nM2 = n - 2;
  43.   pt1.x = nodes[0].x;
  44.   pt1.y = nodes[0].y;
  45.  
  46.   for (i = 1; i < nM2; i++) {
  47.  
  48.     /* B-spline: cubic approximating spline */
  49.     if (splineType == 1) {
  50.       for (j = 0, t = 0.0; j < nT; j++, t += deltaT) {
  51.         t1 = -t * t * t + 3.0 * t * t - 3.0 * t + 1.0;
  52.         t2 = 3.0 * t * t * t - 6.0 * t * t + 4.0;
  53.         t3 = -3.0 * t * t * t + 3.0 * t * t + 3.0 * t + 1.0;
  54.         t4 = t * t * t;
  55.         div = 6.0;
  56.         xOut = (t1 * (double) nodes[i - 1].x + t2 * (double) nodes[i].x
  57.               + t3 * (double) nodes[i + 1].x + t4 * (double) nodes[i + 2].x)
  58.           / div;
  59.         yOut = (t1 * (double) nodes[i - 1].y + t2 * (double) nodes[i].y
  60.               + t3 * (double) nodes[i + 1].y + t4 * (double) nodes[i + 2].y)
  61.           / div;
  62.         pt2.x = (long) (xOut + 0.5);
  63.         pt2.y = (long) (yOut + 0.5);
  64.         if (INIMAGE (pt1, pt2))
  65.           drawline8 (image, imgSize, pt1, pt2, 0);
  66.         pt1.x = pt2.x;
  67.         pt1.y = pt2.y;
  68.       }
  69.     }
  70.  
  71. /* cardinal spline: cubic interp. spline */
  72.     else {
  73.       for (j = 0, t = 0.0; j < nT; j++, t += deltaT) {
  74.         t1 = -t * t * t + 2.0 * t * t - t;
  75.         t2 = 3.0 * t * t * t - 5.0 * t * t + 2.0;
  76.         t3 = -3.0 * t * t * t + 4.0 * t * t + t;
  77.         t4 = t * t * t - t * t;
  78.         div = 2.0;
  79.         xOut = (t1 * (double) nodes[i - 1].x + t2 * (double) nodes[i].x
  80.               + t3 * (double) nodes[i + 1].x + t4 * (double) nodes[i + 2].x)
  81.           / div;
  82.         yOut = (t1 * nodes[i - 1].y + t2 * nodes[i].y
  83.                 + t3 * nodes[i + 1].y + t4 * nodes[i + 2].y) / div;
  84.         pt2.x = (long) (xOut + 0.5);
  85.         pt2.y = (long) (yOut + 0.5);
  86.         if (INIMAGE (pt1, pt2))
  87.           drawline8 (image, imgSize, pt1, pt2, 0);
  88.         pt1.x = pt2.x;
  89.         pt1.y = pt2.y;
  90.       }
  91.     }
  92.   }
  93.  
  94.   return (0);
  95. }
  96.